/** Author: David Powell
This file produces the figures for Appendix F
**/


clear all
set more off
set mat 800
set seed 8721

global dir "/jules/b/dpowell"
global DATA "${dir}/purdue/replication/DATA"
global OUTPUT "${dir}/purdue/replication/output"


use ${DATA}/finaldata

gen overdose_rate=100000*(overdose/totpop)
gen overdose_45t54_rate=100000*(overdose_45t54/totpop)
gen death_rate=100000*deaths_nhw_45t54/totpop

tab stf, gen(ss) nof
tab year, gen(ttt) nof
drop ttt13
foreach nnn of numlist 1983/1994 1996/2017 {
	gen tt_`nnn'=nontripl*(year==`nnn')
}


gen beta1=0
gen low1=0
gen high1=0
gen pred=death_rate
qui reg overdose_45t54_rate ss* tt_* ttt*  [aw=totpop] if year>1989, cluster(stf) 
qui foreach nnn2 of numlist 1990/1994 1996/2017 {
	replace beta1=_b[tt_`nnn2'] if year==`nnn2'
	replace pred=death_rate-_b[tt_`nnn2']*(114/(114-43)) if nontripl==1 & year==`nnn2'
	replace pred=death_rate-_b[tt_`nnn2']*(43/(114-43)) if nontripl==0 & year==`nnn2'
	boottest tt_`nnn2', boottype(wild) weighttype(webb) reps(9999) seed(23857389)
	matrix A=r(CI)
	replace low1=A[1,1] if year==`nnn2'
	replace high1=A[1,2] if year==`nnn2'	
}

gen counter=overdose_rate
qui reg overdose_rate ss* tt_* ttt*  [aw=totpop], cluster(stf) 
qui foreach nnn2 of numlist 1996/2017 {
	replace counter=overdose_rate-_b[tt_`nnn2']*(114/(114-43)) if nontripl==1 & year==`nnn2'
	replace counter=overdose_rate-_b[tt_`nnn2']*(43/(114-43)) if nontripl==0 & year==`nnn2'
}

preserve
gcollapse death_rate pred overdose_rate counter [aw=totpop], by(year) fast
replace pred=death_rate if year<1996


# delimit ;
twoway (connected death_rate year if year>1989 , color(black) msize(vtiny) lwidth(thick) sort yaxis(1))   (connected pred year if year>1994 , color(green) lpattern(dash) msize(vtiny) lwidth(thick)  sort yaxis(1)), graphregion(color(white)) xline(1996) xlabel(1990(2)2017, angle(45))  legend(off) 
scheme(s2mono)  ylabel(,nogrid) legend(off)
xtitle("") ytitle("All-Cause Mortality per 100,000" ,  axis(1)) 
text(415 2003 "Mortality Rate", place(e))
text(392.5 2003.7 "Counterfactual Mortality Rate", place(e))
;
gr export figF2b.eps, replace;


twoway  (connected counter year if year>1994 , color(green) lpattern(dash) msize(vtiny) lwidth(thick)  sort yaxis(1))
(connected overdose_rate year if year>1982 , color(black) msize(vtiny) lwidth(thick) sort yaxis(1))  , graphregion(color(white)) xline(1996) xlabel(1983(2)2017, angle(45))  legend(off)
scheme(s2mono)  ylabel(,nogrid)
xtitle("") ytitle("Overdoses per 100,000" ,  axis(1))
text(17 2005.2 "Drug Overdose Rate", place(e))
text(2.6 2001.2 "Counterfactual Drug Overdose Rate", place(e))
;
gr export figF1.eps, replace;



# delimit cr
restore

preserve
bys year: keep if _n==1
# delimit ;
twoway(rarea low1 high1 year if nontripl==1,  sort color(gs14))  (connected beta1 year, color(black) msize(vtiny) sort yaxis(1)) 
, graphregion(color(white)) xline(1996) yline(0, lpattern(dash) lcolor(gs6)) xlabel(1990(2)2017, angle(45))  legend(label(2 "Estimate") label(1 "95% Confidence Interval")  rows(1) order(2 1))  ylabel(,nogrid)
scheme(s2mono) xtitle("") ytitle("Coefficient Estimate" ,  axis(1)) ;
gr export figF2a.eps, replace;
# delimit cr
restore


